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ABSTRACT 

The probability that a particle, crossing the shock along a given direction, 
be reflected backwards along another direction, was shown to be the key ele¬ 
ment in determining the spectrum of non-thermal particles accelerated via the 
Fermi mechanism around a plane-parallel shock in the test-particle limit. Here 
an explicit equation for this probability distribution is given, for both the up¬ 
stream and downstream sections. Though analytically intractable, this equation 
is solved numerically, allowing the determination of the spectrum in full general¬ 
ity, without limitation to shock speed or scattering properties. A number of cases 
is then computed, making contact with previous numerical work, in all regimes: 
Newtonian, trans-relativistic, and fully relativistic. 

Subject headings: shock waves - cosmic rays 


1. Introduction 

It is possible to formulate particle acceleration around shocks, in the test particle limit, 
in an exact form, for arbitrary shock speeds, and arbitrary anisotropies in the particle dis¬ 
tribution function (Paper I, Vietri 2003). It turns out that the physically relevant boundary 
conditions are those at the shock, not those at downstream infinity. These conditions at the 
shock fix the particle spectrum which was shown to be (under some simplifying assumptions) 
an exact power-law in the particle momentum, independent of the shock speed and of all 
assumptions about the scattering properties of the medium. 

The exact value of the spectral index depends in a detailed way upon the diffusion prop¬ 
erties of the media involved. It was shown in Paper I that this dependence is encapsulated 
in two functions, P„(/ro,/^) and Pd(/io,/i) respectively, which denote the conditional prob¬ 
abilities that a particle entering the upstream (downstream) medium along a direction /io. 
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will leave it for the downstream (upstream) section along a direction /i^. In paper I (but see 
also Freiling, Vietri and Yurko 2003 for some tricky mathematical details), these functions 
were built from the eigenfunctions of the angular part of the scattering problem. Though 
conceptually satisfactory, this procedure is however computationally tricky, and one may 
wonder whether, given the scattering law for the particles, it may be possible to determine 
an equation for the scattering probabilities that yields them directly, without having to build 
the infinitely many eigenfunctions of the angular equation. It turns out that this is actually 
possible, by drawing upon the analogy with scattering atmospheres (Chandrasekhar 1949), 
and by using an obvious symmetry property of the problem at hand. 

It is the purpose of this paper to present this derivation, and a handful of illustrative 
computations, with the major aim of showing the ability of the method to reproduce well- 
established results, drawn from the literature. A more detailed analysis of the spectra in 
several interesting cases, including the important hyperrelativistic regime, will be presented 
elsewhere. 

The plan of the paper is as follows: in section 2, we present the derivation of the 
equations which determine Prf(/io,/i) and Pt,(/io,/i)- In section 3, we show that a diffusive 
equation identical to that originally given by Kirk and Schneider (1987) is recovered in the 
Fokker-Planck regime (or small-pitch-angle scattering, SPAS) of our scattering equation, 
Eq. 2 . In section 4, we describe our full method to compute numerically the spectral slope 
and angular distribution function. In section 5 we show that our method reproduces very 
well results which have already appeared in the literature, in the newtonian, trans-relativistic 
and fully relativistic regimes (the ultra-relativistic regime will be discussed in forthcoming 
publication). We conclude in section 6. 


2. The equations for P^ and 

In paper I, we showed that a suitable, relativistically covariant (though not manifestly 
so) equation for the transport of the particle distribution function for arbitrary pitch angle 
scattering in the presence of a large-coherence length magnetic held is: 

Off Of 

7(m + /X)^ = y i-Weifi', /i, 0, 0) + We(/i, /i', 0, 0')/(/^', 0')) 4^'#' + • (1) 

In this equation, is the distance from the shock along the shock normal, in the shock frame. 
All other quantities are measured in the huid frame: u is the huid speed with respect to 


^In these expressions, all angles are measured in the fluid frame. 



3 


the shock, in units of c, and 7 = \I\J\ — v? its associated Lorentz factor, / is the DF of 
particles with impulse p in the fluid frame (which, by assumption, does not change because 
we assumed pitch-angle scattering), and We{p, is the scattering probability that a 

particle, originally moving along a direction making an angle with the shock normal such 
that its cosine is /i', and along 0', be deflected along a new direction p, 0. The fluid is 
assumed to move in the positive 2 ; direction, so that its speed is at all places > 0 , and the 
angles are measured from the direction of motion of the fluid. 

The scattering is assumed linear (thus these results apply only in the test particle 
approximation), but the scattering angle is nowhere assumed small; we shall derive results 
for the diffusive approximation in a later section. This equation is the plane-geometry, fully 
relativistic version of an equation originally developed by Gleeson and Axford (1967) for 
spherically symmetric, non-relativistic shocks. 

Compared with Eq. 7 of Paper I, the injection has been dropped because we are con¬ 
sidering particles at very large energies, much larger than the injection energy; we now also 
specialize to the case in which no long-coherence length magnetic held is present. The treat¬ 
ment of the problem in full generality will be presented elsewhere (Morlino, Blasi and Vietri, 
in preparation). 

Besides dropping the term uj{df /d(j)), we remark that the system has now become 
symmetric for rotations around the shock normal, so that the DF / will be independent 
of the angle 0. Furthermore, the elemental scattering probability We can depend only on 
the angle 0 — 0 ^ between the incoming and outgoing direction, We = We{p,p',(l) — 00 , even 
though it can depend separately upon p, and p'. In this symmetric case, we can average the 
previous equation over either 0 or 0 ', obtaining 

7(M + h)^ = / i-w{p',p)f{p) + w{p,p')f{p'))dp'dl)', (2) 

where 

= ^ f ^ ~ ’ 

an equation which will be used later on. 

2.1. The downstream case 

In order to hx the ideas, we consider the problem of determining Pd, the conditional 
probability for the downstream frame; we shall explain later how to adapt these results to the 
upstream frame. Take the origin of coordinates to be located at the shock, the downstream 


( 3 ) 
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section to he at z > 0 , and define u as the modulus of the shock speed with respect to the 
fluid. We are using variables in the fluid frame: n is thus the cosine of the angle between 
the particle speed and the shock normal, in the fluid frame. This means that the particles 
which can cross from upstream to downstream [i.e., those with a positive component along 
the shock normal of the velocity with respect to the shock) are those with 

u + n > 0 , (4) 

while those that cross the shock backwards, from the downstream to the upstream section, 
have 

u + fi < 0 . (5) 

Now let / represent a beam of particles, all entering the downstream section along the 
same direction, /io; we showed in paper I that the particle flux across a surface hxed in the 
shock frame, but expressed in downstream coordinates, is given by: 

F = Foj = j 7(m + n)fdia . (6) 

We consider now a pencil beam of particles, all moving initially in the direction fio- Thus, 
at the shock: 

/(h) = —- ho) u + n>0, (7) 

u + ja 

with Fo an obvious normalization. It was shown in paper I that Eq. 7 provides a suitable 
boundary condition for Eq. 2, which then returns the outgoing flux; by dehnition, this is 
precisely the required P^. We thus have 

-(Mrf + /r)/(/r) = PoPd(ho,h) u + fi<0. (8) 

To begin, let us dehne the inward (/+ = f{fi)\u+fi>o) and outward (/_ = /(/i)|t,+^<o) 
parts of the distribution function, and let us consider a surface at a distance from the 
shock, hxed in the shock frame. The hux of particles moving backwards toward the shock, 
through this surface, is: 

7(m + /i)/-(/i) = - / /+(/i')7(w + h')^rf(h', h)c^' • (9) 

J —U 

The minus sign accounts for the fact that m + /r < 0, while all terms inside the integral are 
positive. 

Strictly speaking, in the above equation, we should have written Pd{p,',p,,z) instead of 
simply Pd{n',fi). In other words, for a general situation, we cannot assume that the condi¬ 
tional probability Pd be independent of its location inside the downstream region. However, 
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in this case this assumption is fully justihed, because the downstream region is semi-inhnite: 
in other words, it remains identical to itself whenever we add or subtract hnite regions, and 
its diffusive properties also remain unaltered whenever we add or subtract hnite regions. It 
follows that the conditional probability Pd from some hnite to downstream inhnity must 
be exactly identical to that from 0 to downstream inhnity. This invariance principle, and 
this whole discussion, are identical to those in Chandrasekhar (1949, Chapter 4, Sections 
28-29), and provide an exact justihcation for the equation above. 

We now diherentiate the above with respect to 

7(M + /i)^ = -/ + ( 10 ) 

while Eq. 2 gives us 

1 p — u 

u J-1 

l{u +=-d{fi)f+ +j dp,'w{n,fi')f+{n') + j , ( 12 ) 

where we used the shorthand 

d{jj,) = j w{jj,', jPjdjj,' . (13) 

We now evaluate Eq. 11 at 2 : = 0, using Eqs. 7 and 8, to obtain 


df- 

7(m + + 


7 (m + /i) 



Fod{fi) 

7(M + h) 


Pdi^fJ'Oj j-d) p 


7 (m -h Ho) 



dpi' 


w{p„H')FoPd{Ho,P') 
7(m + h') 


(14) 


We can now proceed analogously for Eq. 12: we again evaluate it at = 0, and use Eqs. 7 
and 8, to obtain: 


liu + h) 



Fodin) Xi ^ , Fow{h,Ho) f ^ ,w{H,H')FoPd{Ho,H') ^ 

tTW Mu + Mo) * WTT) ■ 


As our last step, we plug Eqs. 14 and 15 into Eq. 10 to obtain: 


/^) 


f diyfJjo) 

+ /^O 


d{H) ^ 

u + pij 


w{h, Ho) ^ f Pdin : 1 Po) 

U p Ho J—u U p Ho 



dn' 


Pd{Ho,P')w{H,P') 
u P h' 



dn'PdiH^P) 



dn" 


w{H',H")Pd{,Ho^p") 
u P h" 


( 16 ) 


which is the sought-after result: an equation relating Pd only to the dihusion probability w. 
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2.2. The upstream case 


A completely analogous derivation holds in this case. The pencil beam is now entering 
the upstream section, so that, at the shock, we have: 

/(h) =-y^<^(h-ho) M + /i<0 (17) 

u + fi 

and we have the obvious dehnition, from paper I; 

/(h) = T’n(ho,h) M + /i>0. (18) 

u + jJ, 

At a distance from the shock the flux of particles moving toward the shock is 


{ j- 


/■““ df- 

j 4i'7(M + p')^«(h',h)^ 


(19) 


where the same argument about the use of P„(/i',p) instead of applies as above. 

We may evaluate the above equation at 2 : = 0 using Eqs. 11 and 12 (which still hold), and 
Eqs. 17 and 18, and plug the result into Eq. 19, to obtain 


\U + fio u + fi J 
,w(/i',/io)Pu(h',h) I 


w(h, ho) _ f ^/'U^(h) h )'^M(ho) h ) _|_ 


dfl 


'-1 


U + /io 


U + fio 
djj,'Puin', fJ,) I dpt,' 


'-1 


u + pt' 
, w{pL',pL")Pu{pio,pi" 
U + pi" 


( 20 ) 


which is our other hnal result. 


3. The small—pitch—angle limit 


In this section, we wish to derive the small pitch-angle scattering (SPAS) limit of our 
scattering equation, which holds for arbitrary deflection angles, not just small ones. Before 
proceeding with the derivation of the Fokker-Planck equation, we need to establish a very 
useful result. The function w is not completely arbitrary, but is subject to the following 
requirement. We expect that, when the distribution function (DF) is isotropic (he., when 
/(/i) = const.), there will be no further evolution: in other words, we expect in this case 
the RHS of the above equation to vanish. For this to occur, we shall not take the obvious 
condition 



w(/i, pi')dpi = 


f+i 




w{pi\ pdjdpt 


( 21 ) 
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but its stronger version, 

w{^, = w{^', H) . ( 22 ) 

The above equation is a statement of detailed balance, which we may assume to hold because 
these scattering processes are exactly those which lead to isotropization of the DF, if thermal 
equilibrium were to hold. In our case we have neglected (for excellent, well-known reasons) 
the exchange of energy of the particles with the background fluid, so we cannot expect true 
thermal equilibrium to hold, but since we have included all major scattering processes, we 
expect equilibrium (he., detailed balance except for energy equipartition) to hold at least in 
so far as the angular part is concerned. It is probably worthwhile to remark that a similar 
assumption on w is made by Chandrasekhar (1949, Ch. IV, Section 31, Eq. 29-1). 

We now rewrite Eq. 2 by means of a new scattering probability, a), to be defined 

as follows: 

W{fj,';a)(h = + a, + a) . (23) 

W is simply the probability that a particle is deflected from its initial direction along /i', to 
a new direction jj,' + a. By means of this new quantity Eq. 2 can be recast as: 

7(w + /^)^ = / db(IF(/i-a;a)/(/i-a)-IF(/i;-a)/(/i)) (24) 

while the detailed balance equation becomes: 

hF(/i — a; a) = hF(/r; —a) . (25) 


The scattering equation in the form 24 is exactly identical to that given by Landau 
and Lifshtiz (1984, Section 21, Eq. 21.1), so that their derivation of the Fokker-Planck 
immediately applies. Please note however that our sign convention differs from theirs: what 
we call a is for them —q. Following step by step their derivation we hnd that 


7(m 


dz ^ dfx 


Af + B^ 

OjJ, 


(26) 


with the dehnitions 


r*H-oo 


A = 


da aW(ju; a) 


dB 

(9/i 


^+00 

B = - da a'^W (/r; a) . 

2 J -OD 


(27) 

(28) 


In the two integrals above, we have taken the integration range to extend from —cx) to -|-cx), 
because we have assumed the validity of the SPAS regime, in which case W{fj,; a) has one 
very strong, narrow peak around a = 0: 



We show now that, as a consequence of detailed balance, A = 0. We hrst notice that 


where /i' = /i + a, because of detailed balance, so that 

|H/(«a) = AwV;a). 


( 29 ) 


(30) 


We use the above into 

dB 1 /■+“ 


d 1 d 

da a‘^ — W(jj,;a) = - da a^—W] a) 
dfi 2 dfi 2 da 


(31) 


which can now be integrated by parts to yield 


dB 

dfi 


r*+oo 


da aWa) 


(32) 


We now use again detailed balance, and the variable change x = —a in the dehnite integral 
above, to obtain 


dB 

dfi 


^+oo 


da aW{fi] a) 


(33) 


which, when inserted into Eq. 27, yields A = 0. We thus obtain, for the scattering equation 
in the SPAS, or Fokker-Planck, limit; 


7 (m + fi) 


dl 

dz 


A 

dfi 



(34) 


with B given by Eq. 28. The use of detailed balance implies that the scattering integral 
reduces, in the SPAS limit, to the divergence of a vector proportional to the gradient of the 
DF, exactly like in the classical heat conduction problem. 


As a simple test of this treatment, we now re-derive a well-known result, he., that, in 
the isotropic limit, 

B = K{l-fi^) (35) 

with K a constant. The isotropic limit means that the scattering coefficient probability W 
can only depend on the angle 9 — 9' between the directions of motion. We shall take thus 


U-(wa) = A'/(0 


(36) 


where 


sm{9 - 9') 
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with a a parameter which we take a 1, in order to enforce the SPAS regime. Also, we 
assume that f{x) is an even function of x, with a single minimum around a; = 0 (again, 
SPAS regime); we also assume, for the moment, that W is normalized to unity: 



daN f[a) 


1 , 


(38) 


and will come back to this assumption in a moment. 

We now compute B from Eq. 28. Since we shall be taking the limit a ^ 0, only the 
argument sin{6 — 6^') ~ cr 1 will contribute to the integral, so that, dehning e = 6 — 6', 
we hnd 

e ^ - (39) 

a 

a = jj! — 11 = cos o' — cos 9 = cos{9 + e) — cos 6* ~ — sin Oe (40) 

so that the normalization condition, Eq. 38 requires 

/*+oo ^ 


= ( sin 9(t 


a a 


(41) 


B then becomes 


B = - 
2 


(* + 00 T ^ + oo 

da a'^W{fi; a) ^-Nsin^ 9 / = 


. 

sm 9(7 —x:—- 

f(x)dc 


(42) 


which is clearly identical to Eq. 35. If we now assume W to have arbitrary normalization 
(because it represents the probability of scattering per unit time), we see that the above 
derivation is changed only through an inessential multiplicative constant. This completes 
our derivation. 


3.1. Scattering probability 

Later on, we shall compare numerical results based upon our new approach with those 
which previous authors have obtained in the fully isotropic, small pitch angle scattering limit. 
In order to do this, we have to determine a suitable form for w in this case. It is clear that 
any function very strongly peaked around the forward direction simulates the SPAS regime, 
like, for instance. 


w (X e 


a <^1 . 


(43) 
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However, not all such peaked functions reproduce the isotropic case, which is dehned as that 
in which the elemental scattering probability We is independent of anything except for the 
deflection angle: 

WeifJ,, /i^ 0, 00 = w{cOS 0) (44) 


with 0 given by 

cos 0 = /i/i' + a/1 — /i^ a/1 — /i'^ cos(0 — 00 • (45) 

We have seen above that the transport equation in the SPAS regime does not depend upon 
the exact form of the scattering function We, so that we are free to choose any analytically 
tractable one. We choose: 


We (cos 0) 



1 —cos © 
a 


’) 


a <^1 , 


(46) 


which is obviously suitably normalized to unity when the integral is made over all scattering 
angles. Using Eq. 3 we hnd now 


1 




Wed{(j) - 00 


1 

-e - /o(-- 

a a 


(47) 


where Io{x) is Bessel’s modihed function of order 0, and use has been made of the formula 
(Gradshteyn and Ryzhik, 1994, 3.339): 


exp( 2 :cosa;)cfc = . (48) 

In the following we shall use Eq. 47 when computing spectra, taking a 1 for the isotropic 
SPAS limit and a 3> 1 for the isotropic large angle scattering regime. 



4. Calculation of the spectrum and angular distribution function 

In this section we outline the procedure to follow in order to calculate the spectral slope 
and the angular distribution function of the accelerated particles, for an arbitrary choice of 
the shock velocity and of the scattering properties of the fluid. 

For the purpose of making contact with previous literature on the subject, we carry 
out our calculations with a scattering function w(/i, jj!) as given in Eq. 47, choosing a small 
value of the parameter a for the case of small pitch angle scattering, and a large value of 
(j to reproduce the regime of large angle scattering. Eq. 47 ensures that the scattering 
remains isotropic in both regimes. However, we wish to stress that we elected to use an 
isotropic scattering in order to make contact with all other authors (so far), who have only 
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used isotropic scattering, even though this is by no means required by our method, which 
can effectively deal with anisotropic scattering. 

It is worth stressing that our method does not require any intrinsic approximations or 
assumptions on the type of scattering; the only approximations are due to round-off in the 
numerical solution of the integral equations for Pd and P„. 

The procedure for the calculation of the slope of the spectrum of accelerated particles 
is as follows: for a given Lorentz factor of the shock ( 7 s/i), the velocity of the upstream 
fluid u = Psh is calculated. We then determine the velocity of the downstream fluid Ud from 
the conventional jump conditions at the shock, after assuming some equation of state for 
the post-shock fluid. Given u, Ud and a scattering function w, we can solve the integral 
equations for Pu and Pd iteratively. We recall here that these are functions of an entrance 
and an exit angle, and the convergence of the method can be numerically time consuming, in 
particular for small values of a. Aside from this numerical issue, one should also remember 
that our equations for Pu and Pd are analogous to the H-equation of Chandrasekhar (which 
is however in one dimension only) which is known to admit multiple solutions. It has been 
shown however (Kelley 1982) that iterative methods always converge to the only physically 
meaningful solution, although no information is provided on how fast the convergence is 
achieved. 


An important point to keep in mind is that the probability of return from the upstream 
section is strongly constrained through the condition that all particles must return to the 
shock, being advected with the fluid. It follows that 


d/i P^i(/iQ, jj ,) 


1 


for any value of the entrance angle po- H should be stressed however that this integral con¬ 
straint is not imposed by hand. In other words, though dictated by physical considerations, 
it is automatically implemented by the integral equation for (Eq. 20): it thus provides 
an important numerical check of our computing accuracy. 


The slow convergence of the iterative solution of Eq. 16 has been partly corrected by use 
of the following ruse. The conditional probability Pd obeys the following integral constraint. 
From Eq. 2 we see that at downstream inhnity the condition df/dz = 0 implies that the 
only physically acceptable solution for the distribution function is /(p) = constant = K. At 
downstream inhnity we can therefore write 

f-{ja){u + n) = K{u + ia) = - / d^'U{ia'){u + ia')Pd{ia', ja) = - dia'K{u +ia')Pd{n', n), 

J-Ud J-Ud 

(49) 
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which implies the condition 


du'^Ud + /i) = -{ud + /i) 


' -Ud 


(50) 


for any value of the exit angle /i. Since the function Pd depends only on the scattering 
properties of the fluid, which are assumed to be independent of z, Eq. 50 also applies to any 
including z = 0 which is the position of the shock. 


This condition is of course automatically satisfied at the end of our calculations, while 
it is only approximately satisfied during the iterative solution of Eq. 16. We have found out 
that, by correcting the normalization of Pd so that it obeys the previous integral constraint 
to hold at the end of every iteration, the convergence is much accelerated, above all in the 
SPAS regime. We are unable to explain why this occurs, but shortening of the computation 
times by a factor of 3, or more, have thus been achieved. 


Once the two functions Pu and Pd have been calculated, the slope of the spectrum, as 
discussed in paper I, is given by the solution of the equation: 


{ud + n)g{n)= [ d^Q'^{^,n){ud + 09{0, ( 51 ) 

d-Ud 

where 

[ diyPu{iy,fi)Pd{^,iy)(^ —_ ^'52) 

J-l yl-Urel^J 

Here Urei = is the relative velocity between the upstream and downstream fluids 

and g{jj) is the angular part of the distribution function of the accelerated particles, which 
contains all the information about the anisotropy. A point to notice concerns the function 
Pu in Eq. (52): all variables and functions here are evaluated in the downstream frame, 
while the Pu calculated though Eq. (20) is in the frame comoving with the (upstream) fluid. 
The Pu that need to be used in Eq. (52) is therefore 


P„(z/,/i) = Pu(i>,/i) 


djl 

djj, 


Pu{,T>,g) 


1 ^rel 
(1 "Uj-eZ/i) 


The solution for the slope s of the spectrum is found by solving Eq. 51. In general, this 
equation has no solution but for one value^ of s. Finding this value provides not only the 
slope of the spectrum but also the angular distribution function gi^jP). 


^The fact that the solution is unique could not be demonstrated rigorously, but in our calculations we 
never found cases of multiple solutions. 
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5. Comparison with previous results 

In this section we apply the method detailed above to several cases that have already 
been discussed in the literature, in order to show that the approach can be successfully 
applied to a variety of cases. More specihcally, we consider the following test situations: 

• Since our method is not based on the assumption of small pitch angle scattering, we 

can demonstrate that the slope of the spectrum of particles accelerated at Newtonian 
shocks is independent of the scattering function although the functions Pd 

and Pu in general do depend on the choice of w. 

• To our knowledge, the only calculation of the functions Pd and Pu which applies to 
shocks of arbitrary velocity for the case of large angle scattering is that of Kato & 
Takahara (2001). We compare our results for Pd and Pu with those presented in 
that paper for a scattering function which reproduces the physical case of large angle 
scattering. 

• The slope of the spectrum of accelerated particles in the relativistic regime for the case 
of large angle scattering was compared with the results of Monte Carlo calculations 
carried out by Ellison, Jones and Reynolds (1990). 

• The case of small pitch angle scattering was assumed in most previous calculations. 
We test our approach in this regime in two ways: 1) by showing that the slope of the 
spectrum of accelerated particles in the SPAS approximation, as computed through our 
method for several values of the shock Lorentz factor, compares well with the results 
presented in (Kirk et ah, 2000); 2) by comparing our angular distribution function with 
that obtained by Kirk & Schneider (1987) for a given relativistic shock velocity. 


5.1. The Newtonian case 

The theoretical approach presented here allows us to determine the spectrum and the 
anisotropy in the particle distribution for arbitrary shock velocity and arbitrary scattering 
properties of the fluid. In this section we show that we can reproduce the results expected 
in the non-relativistic limit. We calculated the spectral slope and the angular part of the 
distribution function for several values of the shock velocity and shock compression factor 
r. We also calculated the return probability from the downstream and upstream sections. 
The expected return probability in the Newtonian limit is P^^^ = 1 — 4ud, very accurately 
reproduced by our calculation, for instance for Ud = 10“^ {Pret = 0.96) and for Ud = 2x 10“^ 
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Fig. 1.— Probability functions Pd{fio, n) and P„(/io,Ai) for u = 0.03 and Ud = 0.01 for 
a = 0.01, at fixed values of /xq (as indicated). 


i^ret — 0.92). In all cases that we considered, the condition 
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is conhrmed for any value of /xq within a part in 10^. As a matter of fact, we checked that 
the above equality is satished in all cases we studied, relativistic or otherwise. 

The slope of the spectrum of accelerated particles for n = 3 x 10“^ and Ud = 10“^ 
(compression factor 3) as calculated with our approach is 4.5, in agreement with the theo¬ 
retical prediction (Bell 1978). We obtain the same slope in both cases of small pitch angle 
scattering with a = 0.01 and large angle scattering, that we simnlate by choosing a = 100 
(but any large value gives the same result). The fact that in the newtonian limit the slope of 
the spectrnm does not depend on the scattering properties of the fluid is also in agreement 
with the theoretical prediction of Bell (1978). It is worth stressing that this nniversality in 
the spectrum of the accelerated particles does not reflect in the nniversality of the retnrn 
probabilities Pd and P„. In Fig. 1 we plot Pd^Ho, jj) and Pu{no, IJ') as fnnctions of jj, for the 
values of /io indicated in the hgure and assnming u = 0.03, = 0.01 and a = 0.01. The 

same functions are plotted in Fig. 2 for the case of large angle scattering [a = 100). 


5.2. The probability functions Pd and P^ for the relativistic regime and large 

angle scattering 

As stressed above, onr method is based on the introdnction of the two probability 
distributions Pd and Pu, that can be calculated by solving two integral equations, Eqs. (16) 
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Fig. 2.— Probability functions Pd{fio, fj,) and P„(/io,At) for u = 0.03 and Ud = 0.01 for 
a = 100 (large angle scattering), at fixed values of /xq (as indicated). 



Fig. 3.— Return probabilities Prf(/io,/u) (left) and Pu(/io,/i) (right) that a particle entering 
the downstream (upstream) region with an angle having cosine fio, exits it with an angle 
having cosine p, for the values of po as indicated. 

and (20). The only other approach that we are aware of, that introduces these two functions 
is that presented by Kato & Takahara (2001), based on the theory of random walk. We 
apply our calculations to reproduce Fig. 6 of their paper, obtained for a fluid speed u = 0.5 
both upstream and downstream and in the regime of large angle scattering. 

We plot our results in Fig. 3 for the functions Pd (left panel) and P„ (right panel). We 
carry out our calculations for exactly the same values of po as indicated in Fig. 6 of Kato 
& Takahara (2001) and reported in our Fig. 3 as well. It is worths stressing that for the 
downstream section of the fluid + ho > 0 and Ud + fi < 0, while for the upstream section 
these relations have opposite signs. 
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There is perfect agreement between our results and those of Fig. 6 of Kato and Takahara 
(2001), which is further evidence that our method is perfectly able to reproduce previous 
results, despite the absolute absence of any physical approximation. 


5.3. The spectrum for the case of large angle scattering 

The spectrum of the accelerated particles for the case of pure large angle scattering 
was calculated by Kato & Takahara (2001) and there compared with the previous Monte 
Carlo result of Ellison et ah (1990). The case of mixed small and large angle scattering was 
investigated by Kirk & Schneider (1988). 

We apply our method to determine the spectrum of accelerated particles for pure large 
angle scattering and compare our results with those previously obtained by Kato & Taka¬ 
hara (2001) and by Ellison et ah (1990), which were slightly discrepant for relativistic shock 
speeds. The spectral slopes calculated according with the procedure dehned in the previous 
section are plotted in Fig. 4 for shock speeds u = 0.8 and u = 0.9 as functions of the com¬ 
pression factor r at the shock. The crosses and diamonds are the results of our calculations, 
while the continuous curves represent the interpolation provided by Ellison et ah (1990) to 
their Monte Carlo results. 

It is remarkable that, while the probability distributions Pd and computed by us and 
by Kato and Takahara (2001) agree perfectly, our slopes differ from theirs. We speculate 
that this is due to the different method of computation. Kato and Takahara (2001) use 
Peacock’s recipe: 

P<G >"-3= 1 (53) 

where < G > is the average energy fractional gain when a whole cycle upstream —down¬ 
stream —upstream is performed. However, in paper I, we showed that the correct formula 
is 

P < >= 1 (54) 

which obviously differs from the previous one. The fact that our computations, on the one 
hand, match the functions Pd and Pu of Kato and Takahara, but on the other one reproduce 
so perfectly the slopes of Ellison, Jones and Reynolds (1990) (which use neither our approach 
nor Peacock’s) supports our claim that our formula is the correct one. 
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Fig. 4.— Slope of the spectrum of accelerated particles for a shock with speed u = 0.8 
(diamonds - upper curve) and u = 0.9 (crosses - lower curve) in the regime of large angle 
scattering, as a function of the compression ratio. The continuous lines are the £t to the 
results of Ellison et al. (1990). 


5.4. The SPAS limit 

Most calculations previously presented in the literature adopt the assumption of SPAS, 
in order to reduce the transport equation (Eq. 2) to a Eokker-Planck eqnation, as shown in 
section 3. In the approach presented in Paper I and here, the scattering occurs with arbitrary 
characteristics, embedded in the scattering function w{fi, jj,'). The case of a tc independent of 
pitch angles is that of strong scattering, which is expected to be achieved for large amplitude 
turbulence. The case of small pitch angle scattering in the assumption of isotropic scattering 
can be modeled in our approach by taking a peaked shape for the scattering fnnction, as 
explained in section 4, with a <^1. 
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We expect that in the limit of SPAS the result should not depend upon the detailed 
form of the scattering function, provided this is strongly peaked (see Section 3). 

For a given Lorentz factor of the shock ( 7 ^/ 1 ), the velocity of the upstream fluid u = Psh 
is calculated. For the purpose of comparing our results with those of Kirk et ah (2000), 
we adopt a Synge (1957) equation of state for the downstream fluid and assume that the 
magnetic held is not dynamically important. These are the same assumptions as in (Kirk et 
ah 2000 ). 

We can therefore calculate the velocity of the downstream huid Ud from the conservation 
of mass, momentum and energy at the shock surface. The values of u and Ud are given in 
Table 1. Given these velocities and a width a for the scattering function, we can solve the 
integral equations for and Pd iteratively. 


'^shPsh 

u 

Ud 

slope 

0.04 

0.04 

0.01 

4.00 

0.2 

0.196 

0.049 

3.99 

0.4 

0.371 

0.094 

3.99 

0.6 

0.51 

0.132 

3.98 

1.0 

0.707 

1.191 

4.00 

2.0 

0.894 

0.263 

4.07 

4.0 

0.97 

0.305 

4.12 

5.0 

0.98 

0.311 

4.13 


Following the procedure outlined in section 4, we determine the slope of the power law 
spectrum and the angular part g{fi) of the distribution function. 

The slopes predicted by our calculations for different values of the shock velocity are 
shown in Table 1 for a = 0.01 and agree well with those of Kirk et al. (2000). In Fig. 5 we 
also plot the distribution functions obtained for all the cases reported in Table 1 (the upper 
panel refers to the first four lines in Table 1, while the lower panel refers to the relativistic 
regime, namely the last four lines in Table 1). 

This hgure shows the expected phenomenon of increasing anisotropy in the distribution 
funtion when the shock speed increases: for non-relativistic speeds the distribution function is 
almost perfectly independent of the pitch angle p, but it becomes more and more anisotropic 
in the trans-relativistic regime and in the fully relativistic regime. 

A necessary condition for the SPAS regime to be at work is that a l/dy^^. This 
implies that the SPAS assumption requires increasingly lower values of a when the Lorentz 
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Fig. 5.— Upper panel; DF for ■jshPsh = 0.04 (solid line), 'jshPsh = 0.2 (dotted line), 
IshPsh = 0.4 (dashed line) and 'jshPsh = 0.6 (dash-dotted line). Lower panel; DF for 
IshPsh = 1 (dash-dotted line), 'jshPsh = 2 (dashed line), ■jshPsh = 4 (dotted line) and 
IshPsh = 5 (solid line). 

factor of the shock increases. For a = 0.01, the maximum Lorentz factor for which we can 
assume to be in the SPAS regime is 'jsh ~ 5, and in fact we can already see that the slopes 
that we obtain for 'jshPsh = 5 are slightly smaller than those of Kirk et al. (2000). This 
effect will be discussed in detail in a forthcoming paper. 

It is worth reminding that the SPAS approximation is in fact expected to be broken in 
Nature in the relativistic regime; this is due to the fact that a FP equation describes well 
the scattering only when particles scatter many times within a given angle. In the upstream 
section, particles are caught up by the shock as soon as their deflection angle with respect to 
the normal to the shock is larger than l/'y.sh- For relativistic shocks, this quantity becomes 
small and eventually comparable with the average deflection angle per unit length, and a 
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description by means of a Fokker-Planck equation becomes inadequate, as remarked many 
times in the literature (Kirk et ai, 2000). We stress here that the fact that our method is 
not based upon the SPAS approximation, but on the more general Eq. 2 , makes it more 
widely applicable then previous approaches. Its potential for a more correct investigation of 
the large 7 limit will have to be investigated elsewhere. 

An additional evidence of how important these effects can be for the determination of 
the correct distribution function of the accelerated particles is illustrated below. We calculate 
here the spectrum and the function g{fi) for u = 0.9 and for the relativistic equation of state 
for the downstream gas {uud = 1/3), as adopted by Kirk & Schneider (1987). We carry 
out our calculations for a = 0.05, 0.03, 0.01 and 0.005. The slopes in the four cases are 
very close to each other (s = 4.68, 4.69, 4.71, 4.71) but the angular distribution functions 
(normalized to be unity at the peak) appear slightly different. It is remarkable however that 
when the SPAS regime is approached, the DF converges exactly to that found by Kirk & 
Schneider (1987), as illustrated in Fig. 6 , where we plotted g{fi) for the four cases mentioned 
above: compare our results with Kirk and Schneider’s Fig. 6 . The slope that we obtain is 
also in perfect agreement with Kirk & Schneider (1987). 


6. Conclusions 

In Paper I, a new analytic approach to the calculation of the spectrum of particles accel¬ 
erated at a shock surface moving with arbitrary velocity was proposed. In the present paper 
we described the procedure to calculate the probability functions Pd{fio, fi) and Pu{fio, fi) in¬ 
troduced in Paper I. These two functions allow us to calculate the spectrum and the angular 
distribution function of the accelerated particles at the shock, for arbitrary shock velocity 
and for arbitrary scattering function w. 

We tested the approach in several situations, in order to show that all results previously 
appeared in the literature could be successfully reproduced. In particular the following tests 
have been performed: 

• we conhrmed the universality of the slope of the power law spectrum of accelerated 
particles in the case of Newtonian shocks. More specihcally we proved that such slope 
is roughly independent of the scattering properties of the medium, as described by the 
scattering function w. It should be stressed that such universality does not extend 
to the probability functions Pd and P„, which are instead strongly dependent upon 
the scattering properties. Despite the difference between these probability functions, 
their combination results in an angular part of the distribution function of accelerated 
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Fig. 6.— The angular distribution function g{fi) for u = 0.9 and for the relativistic equation 
of state for the downstream gas uud = 1/3. The solid, dotted, dashed and dot-dashed lines 
refer respectively to a = 0.005, 0.01, 0.03 and 0.05. 

particles which is approximately flat, as expected in the Newtonian limit (isotropy). 

• We checked the correctness of the probability functions Pd and Pu in the case of large 
angle scattering by comparing our results with those of Kato & Takahara (2001). We 
emphasize that these two functions contain all the physical information about the 
spectrum and angular distribution of the accelerated particles. 

• We reproduce exceptionally well the results of the Monte Carlo experiments of Ellison 
et ah (1990) on the slope of the spectrum of accelerated particles as a function of the 
shock compression ratio for different shock speeds, in the large angle scattering regime, 
but not Kato and Takahara’s (2001) results. We ascribed this discrepancy to Kato and 
Takahara’s use of the (incorrect) Peacock’s formula for the spectral slope. 
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• We carried out calculations for the case of a scattering function very peaked as a 
function of the difference between the cosines of the entrance and exit angles. We 
adopted & w in the form given by our Eq. 47. The case of small a is expected to 
describe the limit of isotropic small pitch angle scattering, which is adopted in most 
literature on shock acceleration. We showed that in this limit the transport equation 
(Eq. 2) reduces to a Fokker-Planck equation, identical to that used for instance by 
Kirk & Schneider (1987) and Kirk et al. (2000). We carried out computations for 
a = 0.05 — 0.005, and we compared our predicted spectral slopes with those of Kirk 
et al. (2000) for different values of the shock velocity. Our results are illustrated in 
Table 1 and show that our method reproduces faily well the SPAS results. Moreover, 
our approach shows the expected deviations from the SPAS regime when the condition 
(j 1/(dq^^j) is violated. 

• We have also compared onr distribution fnnction at the shock with that published by 
Kirk and Schneider (1987), obtaining very good agreement, and an identical slope. We 
also show that the distribution fnnction becomes increasingly different when the strict 
SPAS regime is broken (Fig. 6). 

On the basis of the tests carried ont to check onr approach, we can claim that this method 
of calculation of the spectrum and angular distribution fnnction of particles accelerated 
at shocks of arbitrary velocity and arbitrary (but spatially homogeneous) scattering 
properties of the fluid is fully successful in reproducing old results, and it can now be applied 
to situations which go beyond the well stndied case of small pitch angle scattering. In a 
forthcoming paper we will describe the potential of this new approach to describe particle 
acceleration at ultra-relativistic shocks, including the effect of a possible large scale magnetic 
held. 


A sincere thanks is dne to the referee, whose selhess work has considerably improved an 
early version of this work. 
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